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Abstract 

The  object  of  this  study  is  to  derive  a  set  of  equations  which 
predict  the  results  of  orbital  maneuvers  of  vehicles  using  constant  low 
thrust*  These  equations  are  developed  by  simplifying  the  problem  to 
circular  orbits  in  a  two  -  body  dynamic  system.  The  results  are  presented 
in  three  parts.  The  first  part  solves  for  the  equation  of  the  coplanar 
radius  change  problem.  The  second  finds  the  equation  for  the  minimum  fuel 
inclination  change  problem.  The  third  part  of  this  study  puts 'the 
equations  from  the  first  two  parts  together  to  solve  a  minimum  fuel 
transfer  problem  involving  both  radius  and  inclination  changes. 

The  solution  of  the  minimum  fuel  transfer  problem  is  not  to  perform 
the  total  radius  change  and  then  perform  the  total  inclination  change. 
Instead,  the  solution  is  to  perform  both  radius  and  inclination  changes 
on  a  per  orbit  basis. 
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ANALYTICAL  EQUATIONS  FOR  ORBITAL  TRANSFER  MANEUVERS 
OF  A  VEHICLE  USING  CONSTANT  LOW  THRUST 

I  Introduction 

The  launch  of  the  Space  Shuttle  Columbia  ushered  in  a  new  era  in  the 
exploration  and  utilization  of  space.  Because  of  the  reusability  of  this 
vehicle,  the  frequency  and  length  of  spaceflights  will  increase  signifi¬ 
cantly.  As  more  and  more  work  and  research  is  performed  in  Earth  orbit, 
the  need  for  new  and  more  fuel  efficient  propulsion  systems  will  arise. 
Since  February  1970,  one  of  these  new  propulsion  systems  has  been  orbiting 
the  Earth.  This  system,  a  mercury  bombardment  ion  thruster,  is  onboard 
the  SERT  II  (Space  ELectric  Rocket  Test  II)  satellite.  It  was  this 
satellite  which  proved  that  a  constant  low  thrust  could  alter  the  orbit 
of  a  satellite  (Ref  1).  In  the  years  since  1970,  major  advances  in  fuel 
efficiency  and  thruster  service  life  have  been  made  for  electric 
propulsion  systems  (Ref  1,2). 

But  how  can  a  low  level  of  thrust  change  the  orbit  of  a  satellite? 

The  answer  lies  in  the  fact  that  the  thrust  is  constant  over  a  long 
period  of  time.  The  velocity  of  a  four  metric  ton  satellite  thrusting  at 
0.8  Newtons  for  30  days  can  be  changed  by  as  much  as  1,866.0  kilometers 
per  hour.  This  velocity  change  can  be  achieved  with  substantial  fuel 
savings  using  electric  propulsion  instead  of  a  more  conventional  chemical 
propulsion  system- (Ref  1). 

Changing  the  velocity  of  a  satellite  such  as  described  above  can  be 
used  to  boost  a  satellite  into  a  planetary  mission  orbit.  This  is  only 
one  of  the  many  uses  of  electric  propulsion  systems.  Some  of  the  more 
immediate  uses  will  be  as  stationkeeping,  attitude,  and  repositioning  jets 
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for  satellites.  Future  uses  include  the  placement  and  retrival  of  sat¬ 
ellites  in  geosynchronous  orbit  and  the  movement  of  large  fragile  space 
structures. 

As  mentioned  above,  one  of  the  more  immediate  uses  of  constant  lov 
thrust  will  be  in  the  repositioning  of  satellites.  This  repositioning 
maneuver  will  require  changes  in  the  radius  and  the  inclination  of  the 
satellite's  orbit.  Reference  3  (chapter  3)  derives  a  set  of  equations 
which  give  the  change  in  velocity  needed  to  achieve  certain  changes  in 
radius  and  inclination.  But  these  equations  can  not  be  used  for  a 
satellite  using  constant  low  thrust  because  they  were  derived  with  the 
assumption  that  the  vehicle's  velocity  vector  could  be  changed  instan¬ 
taneously.  This  is  a  valid  assumption  for  high  thrust  chemical  propulsion, 
but  not  for  low  thrust  electric  propulsion.  Therefore,  a  set  of  equations 
must  be  derived  for  the  constant  low  thrust  case. 

This  paper  has  three  objectives.  The  first  objective  is  to  derive  an 
equation  which  will  give  the  amount  of  radius  change  per  orbit  of  a  vehicle 
using  constant  low  thrust.  It  will  be  assumed  that  the  vehicle  stays  in 
the  plane  of  its  starting  circular  orbit.  This  makes  the  problem  similar 
to  the  Holmann  Transfer  solved  in  chapter  3  of  Reference  3.  The  solution 
of  the  constant  low  thrust  maneuver  will  be  presented  in  chapter  2  of  this 
paper. 

The  second  objective  is  to  derive  an  equation  which  will  give  the 
amount  of  inclination  change  of  a  vehicle,  similar  to  the  one  mentioned 
above,  using  a  minimum  fuel  control  profile.  It  will  be  assumed  that  the 
vehicle  will  start  and  finish  in  a  circular  orbit  and  only  change  its 
inclination.  Similar  equations  for  the  instantaneous  high  thrust  case 
can  be  found  in  chapter  3  of  Reference  3,  The  equations  for  the  constant 
low  thrust  case  will  be  derived  in  chapter  3  of  this  paper. 
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The  third  objective  is  to  join  the  equations  derived  in  this  study 


for  radius  and  inclination  changes  so  that  the  total  orbital  transfer 
problem  can  be  solved  optimally.  The  problem  is  to  be  solved  by  finding 
the  amount  of  each  type  of  change  vhich  must  be  made  per  orbit  to  reach 
the  end  conditions  in  a  specified  number  of  orbits. 
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II  The  Coplanar  Radius  Chan  fro  Problem 

The  transfer  between  any  two  coplanar  circular  orbits  is  a  very 
simple  and  useful  maneuver.  Today,  this  maneuver  can  be  accomplished  by 
a  number  of  different  ways.  Among  these,  the  most  economical  is  the 
Holmann  transfer  ellipse  (Ref  3:  163-165).  But  these  methods  are  based 

on  two  assumptions  that  are  not  valid  when  the  vehicle  in  question  is 
using  low  thrust.  The  first  assumption  is  that  the  vehicle  can  change 
its  velocity  vector  instantaneously.  This  assumption  is  violated  since, 
with  constant  low  thrust,  it  may  take  hours  or  days  to  significantly 
change  the  vehicle’s  velocity  vector.  The  second  assumption  is  that  the 
vehicle  can  accelerate  against  the  pull  of  the  earth’s  steep  gravity  well. 
Since  the  acceleration  of  today’s  constant  low  thrust  propulsion  systems 
is  on  the  order  of  10  ^g's,  the  second  assumption  is  no  longer  valid. 
Therefore,  with  these  two  restrictions  in  mind,  a  new  maneuver  and  a  new 
equation  describing  the  results  of  the  maneuver  must  be  sought. 

The  search  for  the  new  maneuver  began  by  looking  at  how  the  radius 
of  a  circular  orbit  can  be  changed  without  violating  the  above  restrictions. 
To  overcome  the  strong  pull  of  gravity,  the  vehicle  must  thrust  in  a 
direction  which  is  along  gravity's  weakest  lines.  This  direction  happens 
to  be  perpendicular  to  the  lines  of  force  of  the  earth's  gravitational 
field.  This  direction  is  parallel  to  the  vehicle's  velocity  vector.  It 
can  be  readily  seen  that  if  the  vehicle  thrusts  in  this  direction  the 
result  is  a  spiral.  This  type  of  maneuver  fulfills  the  objective  of  the 
above  maneuvers  while  at  the  same  time  does  not  violate  the  above  restrio- 
tions.  Therefore,  this  will  be  the  maneuver  used  to  change  the  radius  of 
a  circular  orbit  using  constant  low  thrust. 
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Now  that  the  maneuver  has  been  defined,  an  equation  must  be  derived 
which  will  describe  the  results  of  such  a  maneuver.  If  the  thrust  is  low 
enough,  the  spiral  can  be  modeled  as  a  series  of  concentric  circular 
orbits  with  a  discrete  change  in  the  radius  at  the  point  on  the  circular 
orbit  where  the  thrusting  was  started.  By  looking  at  two  of  the 
concentric  circular  orbits  and  the  curve  that  joins  them  (see  Fig  1),  the 
equation  which  gives  the  amount  of  radius  change  can  be  found. 


Fig  Is  Spiral  Approximation 

The  equation  for  the  specific  energy  of  a  circular  orbit  is 


given  by 


(D 


where  r  is  equal  to  the  radius  of  the  circular  orbit.  The  time  derivative 
of  Eq  1  can  be  gotten  by  using  elemental  differential  calculus.  This 
derivative  is  simply 

d  (EJ  .  d(r)  (2) 

<)t  2rx  at 
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Since  the  energy  of  the  circular  orbit  can  only  be  increased  by  thrusting 
in  the  direction  parallel  to  the  velocity  vector,  the  tine  derivative  of 


E  can  be  written  as 
c 

KE,)  .  F  •  V 

dt  "  «> 

where  F  is  a  specific  force  (force/unit  mass)  and  V  is  a  velocity. 
Equations  1  and  2  can  be  equated  and  the  r  tern  solved  for  to  get  the 
following  equation: 

ill!  =  (F-  v)  U) 

dt  s*  <4) 


The  first  term  of  Eq  4  can  be  evaluated  as  the  magnitude  of  the  two 
vectors  multiplied  by  the  cosine  of  the  angle  between  them.  Since  the 
thrust  vector  is  aligned  with  the  velocity  vector,  the  angle  between 

mm  mm 

them  is  zero  and  the  cosine  term  becomes  one.  Also  the  F  and  V  terms  can 


be  evaluated  as  the  thrust  per  unit  mass  and  the  velocity  of  a  circular 
orbit,  respectively.  This  can  be  done  since  the  distance  between  the 
initial- and  final  circular  orbits  is  very  small.  Inserting  these 
expressions  into  Eq  4  gives  the  following  formula: 


dir)  -  IT 
dt 


(5) 


The  next  step  in  the  solution  of  this  problem  is  to  separate  Eq  5  and 
integrate  each  side  with  the  proper  limits.  This  leads  to  the  solution 
form  which  gives  the  change  in  radius  using  constant  low  thrust  as  a 
function  of  the  thrust  per  unit  mass  and  the  time  set  to  achieve  the 
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change.  This  solution  equation  is 


T 


ivto 


(6) 


Since  the  spiral  has  been  modeled  as  a  circular  orbit  with  a  discrete 
change  in  the  radius  at  the  point  where  the  thrusting  was  started,  r^ 
can  be  written  in  the  following  form: 

-  r,  +  o.r 


This  allows  Eq  6  to  be  written  as 


The  (tp  -  t^)  term  can  be  interpreted  as  the  period  of  a  circular  orbit 
because  the  Ar  is  very  small.  From  Reference  3  (33),  the  period  of  a 
circular  orbit  is 


(9) 


This  allows  Eq  8  to  be  written  in  its  final  form  as 

Ar  ~  T  a  ir  i?*'1  i  l1  ~  no) 

L  *  ftpi 

where  A  r  is  the  radius  change  per  orbit  using  constant  specific  thrust, 
T,  at  radius  r^.  Now  that  the  equation  for  the  coplanar  radius  change 
has  been  found  the  next  step  is  to  find  the  solution  of  the  inclination 
change  problem. 
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Ill  The  Inclination  Change  Problem 

The  second  chapter  of  this  paper  dealt  with  the  transfer  betveen  two 
coplanar  circular  orbits.  But  a  satellite  may  need  more  than  a  radius 
change,  it  may  also  need  an  inclination  change.  It  is  the  object  of  this 
chapter  to  develop  the  equation  which  describes  the  inclination  change. 

Setting  Up  the  Inclination  Change  Problem 

In  chapter  2,  two  restrictions  were  placed  on  the  coplanar  circular 
radius  change  problem.  These  two  restrictions  must  be  observed  whenever 
constant  low  thrust  is  used.  Just  as  in  the  radius  change  problem,  the 
thrust  must  be  oriented  in  a  direction  which  is  along  a  line  perpendicular 
to  the  lines  of  force  of  the  Earth’s  gravity  gradient.  In  the  radius 
change  problem,  this  direction  was  chosen  to  be  parallel  to  the  velocity 
vector  of  the  vehicle.  However,  in  this  problem  the  direction  is  perpen¬ 
dicular  to  both  the  gravity  gradient  and  the  velocity  vector.  Figure  2 
shows  the  orientation  of  the  thrust  vector,  as  described  above,  to  change 
the  orbital  plane  in  a  counterclockwise  direction  about  the  x-axis. 
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The  equation  for  the  inclination  change  problem  can  be  found  by  looking 


at  the  specific  angular  momentum  equations  which  describe  a  two  -  body 
circular  orbit.  Reference  3  ( 16-17)  proves  that  for  a  two  -  body  orbit 
the  specific  angular  momentum  of  the  orbit  is  constant  and  can  be  given  by 

h  n  r  X  V  (11) 

It  also  proves  that  r  and  v  form  a  plane  that  is  fixed  in  space.  This 
plane  is  called  the  orbital  plane.  Since  TS  is  the  vector  cross  product 
of  r  and  v,  h  must  be  perpendicular  to  the  orbital  plane.  The  incline 
tion  of  the  orbital  plane  is  defined  in  Reference  3  (58-59)  as  the  angle,  i, 
between  K  axis  of  the  Geocentric  -  equatorial  coordinate  system  and 


Pig  3s  Geocentric  -  equatorial  Reference  Frame 
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Prom  Fig  4»  it  can  be  seen  that  the  change  in  h  is  defined  as 

Ah  -  -  k, 


(12) 


Fig  4:  Inclination  Change 

The  angle  between  h^  and  h^  is  the  change  in  inclination  between  the  two 
orbit  planes.  As  this  change  in  h  becomes  small,  the  angle  between 
A  h  and  h,  can  be  approximated  by  a  right  angle.  Therefore,  the  angle 
between  h^  and  h^  becomes 


SIN  A I 


Ah 

h 


(13) 


Since  this  equation  holds  only  for  Ah  very  small,  Eq  13  can  be  further 
simplified  by  using  a  small  angle  approximation  for  the  sine  of  Ai» 
Equation  13  can  now  be  written  as 


(14) 
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From  Fig  4,  it  can  be  seen  that  to  change  the  inclination  in  the 
counterclockwise  direction,  a  moment  must  be  applied  about  the  x-axis. 
From  Reference  4  (23),  the  moment  about  the  x-axis  can  be  defined  as  the 
time  derivative  of  the  specific  angular  momentum.  In  equation  form  it  is 


h  =  M 


(15) 


Equation  15  can  be  rewritten  as 

d(K) 

cl*  ' 


(16) 


This  form  can  be  separated  and  substituted  into  Eq  14  to  get  the  following 
form  for  the  change  in  inclination: 

d»  -  (17) 

n 

The  next  step  in  this  derivation  is  to  find  the  equation  for  the  moment 
about  the  x-axis.  From  Fig  2,  it  is  seen  that  the  moment  about  the  x-axis 
is  produced  by  the  component  of  the  thrust  vector  perpendicular  to  the 
orbital  plane.  Therefore,  the  moment  can  be  written  as 


M  =  r  x  T 


(18) 


From  the  notation  in  Fig  2  and  the  definition  of  the  vector  cross  product, 
Eq  18  can  be  written  as 


M  =  vT  si  n  e 


(19) 
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The  substitution  of  Eq  19  into  Eq  17  allows  the  change  in  inclination  to 
be  written  as 


SIN  © 


(20) 


The  change  in  tine,  dt,  of  Eq  20  can  be  written  as 


(21) 


or 


(22) 


The  substitution  of  Eq  22  into  Eq  20  allows  the  change  in  inclination  to 
be  written  as  a  function  of  0  only  since  r  and  v  are  constants  in  a 
circular  orbit.  The  new  equation  for  the  inclination  change  is 


x!t 

KV 


SIN  © 


(23) 


The  next  step  in  this  derivation  is  to  integrate  the  right  hand  side  of 
Eq  23  from  0  to  2 IT  ,  But,  as  stated  above,  r,  h,  and  v  are  all  constar.1 
in  a  circular  orbit.  This  leaves  the  specific  thrust  (thrust/unit  mass) 
in  the  integral. 

It  can  be  readily  seen  that  the  greatest  amount  of  inclination  char, 
is  achieved  by  allowing  the  thrust  to  be  constant  over  the  entire  orbit. 
However,  this  wastes  fuel  because  very  little  inclination  change  is 
achieved  when  the  vehicle  is  close  to  the  x-axis.  If  fuel  usage  is  r.o1 
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a  problem,  then  this  is  the  way  to  achieve  the  greatest  amount  of 
inclination  change  in  the  shortest  amount  of  time.  But  on  most  space- 
flights,  fuel  must  be  carefully  conserved.  Therefore  thrust  histories 
which  are  functions  of  0  are  of  interest.  This  generalization  allows 
Eq  23  to  be  written  as 


.2.TT 

T(e)  sin  © 


(24) 


The  substitution  of  the  definitions  of  h  and  v  from  Reference  3  (29,34) 
allows  Eq  24  to  be  written  in  its  final  form  as 


Q.T 

j  T(B)  si  N  Q 

o 


(25) 


If  the  form  of  T(6)  is  specified,  Eq  25  can  be  integrated  to  find  Ai. 
Since  fuel  has  to  be  conserved,  the  logical  form  for  T(0)  needs  to  be 
chosen  to  minimize  the  fuel  consumption.  This  form  for  T(9)  can  be  found 
by  solving  an  optimal  control  problem  where  the  performance  index,  fuel 
consumption,  to  be  minimized  is  given  by 


The  minimum  of  Sq  26  without  a  constraint  is  obviously  zero.  Therefore 
Eq  25  with  Ai  specified  is  used  as  a  constraint.  This  optimal  control 
problem  is  solved  in  the  following  sections  of  this  chapter. 
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Suboptitnal  Control  Approach 


The  solution  of  all  but  the  simplest  optimal  control  problems  is 
very  difficult  to  find.  Therefore  certain  methods  have  been  developed 
to  approximate  the  solution  of  the  optimal  control  problem  (Ref  5) .  The 
method  used  to  solve  the  inclination  change  problem  is  the  suboptimal 
control  method. 

Suboptimal  Control  Problem.  The  suboptimal  control  problem  is  a 
parameter  optimization  problem.  The  transformation  from  optimal  control 
to  suboptimal  control  is  achieved  by  assuming  that  the  optimal  controls 
of  the  problem  in  question  can  be  approximated  by  a  known  mathematical 
form  having  a  set  of  unknown  parameters.  When  this  transformation  is 
used,  the  performance  index,  J,  becomes 

=  3(A) 


where  A  is  the  set  of  unknown  parameters  associated  with  suboptimal 
approximation.  The  optimal  A*s  are  found  by  solving 

(28) 

Solution  Method.  There  are  many  methods  used  to  solve  the  parameter 
optimization  problem  stated  in  the  preceding  section.  The  method  chosen 
to  solve  the  inclination  problem  is  the  Conjugate  -  Gradient  method  (Ref  6). 

The  Conjugate  -  Gradient  method  is  a  quasi  -  second  order  technique 
for  solving  parameter  optimization  problems.  It  is  very  similar  to  the 
Method  of  Steepest  Descent  otherwise  known  as  the  gradient  method.  The 
only  major  difference  between  the  two  is  that  the  Conjugate  -  Gradient 


J(J) 

1A 


method  carries  along  previously  calculated  gradient  directions  in  a  fading 
memory  form.  It  is  this  fading  memory  that  gives  this  method  its  second 
order  characteristics. 

The  iterative  algorithm  for  this  method  is  as  follows: 

1)  Guess 

2)  Compute  V  J(A^) 

3)  Let  SJ’  =-V  j(Aj) 

4)  Compute  A^+"*  =  A^*  +  c<  ^ 

5)  Compute  V  J(A^) 

6)  If  V  J(A^+^)  —  Tolerance  :  Stop 

7)  Compute  SJ+1  =-Vj(A3+1)  + 

8)  Increment  j,  go  to  step  4 

There  are  several  terms  in  the  above  algorithm  that  need  to  be  defined. 
The  first  is  the  o(  *  term  of  step  4.  This  term  is  a  step  size  that  is 
computed  with  a  1  -  dimensional  search  method  found  in  Reference  7  (55). 

This  search  method  finds  the  minimum  of  a  parabola  in  the  J  vs  <4  plane. 

The  S  terms  are  directions  along  which  the  algorithm  searchs  for  the 
optimal  A’s.  The  f3  ^  term  is  the  fraction  which  tells  the  algorithm  how 
much  of  the  previous  S's  should  be  used.  This  is  the  fading  memory  of  the 
Conjugate  -  Gradient  method.  The  {3  ^  term  is  determined  by  the  following 
formula.  .  2 

ft3  -  )  I  (29) 

P  W3(^)f 

It  should  be  noted  that  the  same  procedure  which  gives  this  technique  its 
second  order  characteristics  can  also  slow  the  convergence  of  the  algorithm. 
Therefore,  the  |S  ^  is  set  to  zero  periodically  to  prevent  the  algorithm  from 
carrying  too  much  past  information.  A  more  detailed  explanation  of  the 
Conjugate  -  Gradient  method  can  be  found  in  Reference  6. 


Solving  the  Constrained  Inclination  Problem 


Control  Variable  For:?,.  It  is  obvious  that  the  direction  of  the 
thrust  vector  must  be  in  opposite  directions  on  either  side  of  the  x-axis 
to  achieve  the  inclination  change  represented  in  Fig  4,  This  means  that 
T(9)  is  an  odd  function.  Also  the  form  of  T(6)  must  be  symmetric  about 
the  90°  point  of  the  orbit.  It  is  at  this  point  that  the  greatest 
moment,  and  the  corresponding  greatest  change  in  inclination,  is  achieved 
by  the  thrust.  Because  of  these  two  qualifications,  a  Fourier  sine  series 
with  the  even  numbered  terras  removed  was  chosen  to  be  the  approximate 
optimal  control  form  (Ref  8:  35-36,  80;  9),  This  form  is 


oo 


T(e)  = 

k-1 


(30) 


Nov;  that  the  form  for  T(0)  has  been  specified,  the  integration  of  Fq  25 
can  be  performed.  The  evaluat  rr  of  the  constraint  equation  results  in 
all  the  terras  except  the  sin  9  term  being  zero.  This  allows  the  first 
coefficient  of  the  Fourier  sine  series  to  be  solved  for.  This  coefficient 

"  A  - 

“l  “  TT  (31) 


Thus  for  a  given  inclination  change,  the  first  coefficient  of  the  Fourier 
sine  series  is  known.  The  remaining  coefficients,  needed  to  give  the 
minimum  fuel  consumption,  are  unknown.  Since  the  first  cofficicnt  is 
known,  it  can  be  divided  out.  This  results  in  T(©)  having  the  following 
form. 


T(e)  - 


a  \ju 

faTr 


XX 


(32) 
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(33) 


where  tt  is  given  by 

TC  -  S\H©  + 


00  / 

H(a. 


s*n  (ek- 1)©) 


where 


(34) 


With  this  form  for  T(6),  Eq  25  no  longer  needs  to  be  used  as  a  constraint. 
Now  only  the  A’s  need  to  be  found  for  minimum  fuel  usage.  But  because  of 
this  new  form  for  T(0),  a  new  constraint  crust  be  used.  The  first  term  of 
Eq  25  can  be  looked  at  as  the  maximum  thrust,  while  the  second  tern  looks 
like  a  throttle  control.  If  no  constraint  is  put  on  the  throttle,  the 
solution  of  the  optimal  control  problem  would  be  an  infinite  spike  at 
the  90°  point  of  the  orbit.  Therefore,  Eq  33  has  the  following  limits: 

O  £  Tt  £  1  (35) 


Numerical  Solution.  The  algorithm  for  the  Conjugate  -  Gradient 
method  requires  the  calculation  of  the  gradient  of  the  performance  index 
with  respect  to  the  A’s.  These  gradients  are  found  using  the  central 
difference  formula  described  in  Reference  10.  This  formula  has  the 
following  form: 


tjffi  _  » S  A»-0  ~  3  (Azk.i  -  S  Aiu-i) 

2.  S  Ajvj.-t 


(36) 


whereJA^  was  calculated  by 

%  -  1  E-0  A 


(37) 
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V. 

To  find  the  perturbed  J’s,  the  integral  of  Eq  26  has  to  be  evaluated. 

Since  the  integrand  of  this  integral  involves  an  absolute  value,  the 

areas  above  and  below  the  xnaxis  of  the  T  vs  9  curve  have  to  be  summed. 

This  means  that  the  roots  of  the  control  equation  have  to  be  found. 

Two  methods  from  Reference  10  (65-71)  are  used  to  find  the  roots. 

The  first  method  is  the  Bisection  method.  This  is  a  brute  force  method 

which  is  used  to  isolate  the  root.  When  the  root  is  isolated,  the 

Secant  method  is  used  to  improve  the  accuracy  of  the  answer.  The 

Secant  method  is  a  Newton  -  Raphson  method  where  the  derivatives  are 

replaced  with  difference  approximations. 

Included  in  the  root  searching  algorithm  is  an  algorithm  which 

finds  the  points  where  the  value  of  T(©)  is  greater  than  T  .  Since 

D&X 

the  thrust  cannot  be  greater  than  T  ,  the  area  above  this  line  cannot 

max’ 

be  included  in  the  calculation  of  the  performance  index.  Therefore,  when 
the  integral  is  evaluated,  the  value  of  the  integrand  is  constant  between 
the  points  where  the  T(6)  curve  lies  above  the  T  line.  The  area  vhere 

ElclX 

the  integrand  is  constant  is  calculated  by  using  the  formula  for  the  area 
of  a  rectangle.  The  areas  on  either  side  of  the  rectangle  are  calculated 
using  a  Runga  -  Kutta  numerical  integrator  with  a  variable  step  size. 
Since  the  areas  on  either  side  of  the  rectangle  are  the  same,  only  one 
of  them  has  to  be  calculated.  The  total  area  is  twice  this  value  plus 
the  area  of  the  rectangle. 

The  algorithm  was  started  with  only  one  A  in  Eq  33»  The  A  was 
changed  until  the  gradient  of  the  performance  index  was  less  than  or 
equal  to  a  set  tolerance.  The  solution  was  then  said  to  be  converged 
and  another  A  was  added.  A’ s  were  to  be  added  to  Eq  33  until  the 
performance  index  failed  to  change  with  the  new  addition.  But  before 
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this  happened,  the  solution  to  the  optimal  control  problem  uas  suggested 

from  the  partial  results  of  the  program. 

Solution  Results.  The  Conjugate  -  Gradient  program  described  above 

was  run  until  nine  A's  had  been  added  to  the  optimal  control  equation. 

The  plot  of  this  equation  can  be  seen  in  Fig  5,  but  the  program  did  not 

use  the  total  equation.  It  used  the  equation  until  it  went  above  the 

T  line  and  then  followed  the  T  line  until  it  drooped  off  to  follow 

max  max 

the  T(0)  equation  again.  This  form  can  be  seen  in  Fig  6.  If  the  small 

wiggles  are  considered  to  be  noise,  then  the  form  which  emerges  is  a 

square  wave  centered  about  the  90°  point  of  the  orbit.  To  test  this  idea, 

a  nine  term  square  wave  approximation  was  used  to  compute  a  performance 

index.  The  value  of  the  performance  index  of  the  square  wave  was 

substantially  lower  than  the  performance  index  of  the  nine  term  Fourier 

sine  series.  Since  it  takes  up  to  forty  terms  to  accurately  describe  a 

square  wave,  the  performance  index  can  be  drastically  lowered  by  using 

the  square  wave  instead  of  the  Fourier  sine  series.  Therefore,  it  uas 

concluded  that  the  solution  of  the  constrained  inclination  change  optimal 

control  problem  is  a  square  wave  with  T  as  its  maximum. 

max 

Bang  -  Coast  -  Banff  Solution 

The  square  wave  solution  of  the  inclination  change  problem  was 
suggested  by  the  results  of  the  suboptimal  control  problem  outlined  in 
the  preceding  sections.  A  plot  of  the  square  wave  can  be  seen  in  Fig  7. 
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Fig  6:  Constrained  Optimum  Control  Profile 
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Fig  7:  Square  Wave 

To  find  the  amount  of  inclination  change  this  control  profile  gives, 
Eq  25  must  be  integrated  between  the  limits  set  forth  in  the  graph  with 


T(9)  =  Tfflax«  Equation  25  can  be  rewritten  as 


■t  _ - 

i'j  _  t*  io>A» 


T,„  a.u, 

I  I  «>»*>  ©  <Je  —  j  ci©  I 

*  x 


(38) 


To  evaluate  Eq  38,  two  integrals  must  be  solved.  Since  thrusting  on 
either  side  of  the  orbit  achieves  a  change  in  the  inclination,  the 
integral  can  be  simplified  to  the  following  form: 


A\  r 


+  o> 


S>N  © 


de 


(39) 


The  evaluation  of  this  integral  and  the  use  of  trigometric  formulas  to 


simplify  the  integrated  form,  leads  to 


A*.  -  _4r  Tn«, 


(40) 


e 
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where  i  is  the  inclination  change  per  orbit  using  constant  specific 


thrust-  T  -at  radius  r  for  w.  the  thrust  time  angle. 

’  max’  ’  ° 

Now  that  the  equations  for  the  radius  and  inclination  changes  have 
been  found,  the  two  separate  equations  must  be  joined  to  form  an  equation 
which  can  be  used  to  describe  the  total  orbital  transfer  problem. 
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IV  The  Combined  Orbital  Transfer  Problem 

In  the  preceding  two  chapters,  the  problens  of  the  coplanar 
circular  radius  change  and  the  inclination  change  were  solved.  Although 
the  general  combined  orbital  transfer  problem  could  not  be  solved,  the 
solution  of  a  specific  example  problem  was  found.  It  is  the  purpose  of 
this  chapter  to  define  and  solve  this  problem. 

Defining  the  Orbital  Transfer  Problem 

The  orbital  transfer  problem  can  bd  stated  as  follows:  Given  the 
initial  and  final  values  of  the  radius  and  inclination  and  the  number  of 
orbits  inwhich  to  achieve  these  changes,  what  amounts  of  radius  end  incli¬ 
nation  changes  must  be  attained  per  orbit  to  match  the  final  values?  To 
simplify  the  problem,  two  assumptions  were  made.  The  first  was  that  the 
vehicle  is  using  full  thrust  in  either  the  radius  change  direction  or  the 
inclination  change  direction,  but  never  dividing  the  thrust  between  the 
two  directions.  The  second  assumption  was  that  the  direction  of  the 
thrust  vector  can  be  changed  instantaneously. 

There  are  two  methods  which  can  be  used  to  achieve  the  final  radius 
and  inclination.  The  first  is  to  thrust  parallel  to  the  velocity  vector, 
changing  the  radius,  until  the  final  radius  is  achieved.  Then  the  thrust 
is  directed  perpendicular  to  the  orbital  plane,  changing  the  inclination, 
until  the  final  inclination  is  achieved.  The  second  method  is  to  change 
both  the  radius  and  the  inclination  in  a  given  orbit.  It  is  the  second 
method  which  is  the  most  fuel  and  time  efficient.  The  inclination  change 
equation,  Eq  40,  is  in  the  ideal  form  to  handle  this  situation.  But  the 
radius  change  equation,  Eq  10,  is  not.  Therefore,  Eq  10  must  be  modified 
to  allow  for  this  situation. 
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Modification  of  the  Radius  Change  Equation 


Equation  10  was  derived  with  the  assumption  that  the  entire  period 
of  the  circular  orbit  would  be  used  for  the  radius  change.  But  for  the 
total  orbital  transfer  problem,  the  time  spent  changing  the  inclination 
must  be  subtracted  from  the  period  of  the  orbit#  Since  the  inclination 
thrust  time  angle,  w,  of  Eq  40  is  in  radians,  the  equation  for  the 
angular  velocity  of  a  circular  orbit  must  be  used  to  make  the  change  to 

seconds.  Also  it  must  be  realized  that  this  change  has  to  be  multiplied 

/ 

by  four  since  there  are  four  such  times  per  orbit#  The  new  equation  for 
the  time  to  change  the  radius  is 


TP  =  1Pt 


(41) 


Equation  10  can  now  be  written  in  the  more  usable  form  of 


Now  that  the  radius  change  equation  has  been  modified,  the  total  orbital- 
transfer  problem  can  be  solved# 


Solution  of  the  Orbital  Transfer  Problem 

As  stated  in  the  first  section  of  this  chapter,  the  amount  of  radius 

and  inclination  change  per  orbit  must  be  found  to  match  the  final  radius 

and  inclination,  r  and  i  respectively,  after  a  fixed  number  of 
’  nom  nom  ’ 

orbits#  This  will  be  done  by  setting  up  a  suboptimal  control  computer 
program  similar  to  the  one  in  chapter  3  of  this  paper. 
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Solution  Method.  The  algorithm  given  in  chapter  3  for  a  suboptimal 
control  Conjugate  -  Gradient  computer  program  will  still  be  used  for  this 
problem.  But  two  changes  have  to  be  made  to  the  complete  program.  The 
first  is  to  change  the  form  of  the  control  equation.  The  following  form 
will  be  used: 


co(0  =  A„  +  (A  R  f 

k*1 


(43) 


where 

R  =C'r  -  Y.)  /  ( W  -  r\) 


(44) 


The  second  change  is  in  the  form  of  the  performance  index.  The  new 
performance  index  is 

1  ** 

"3"  =  '  rt  **•*»''')  (45) 

The  rf  and  if  terms  in  Eq  45  are  the  calculated  final  radius  and 
inclination  achieved  using  Eqs  4 0,  42,  and  43  in  an  iterative  loop  for 
the  specified  number  of  orbits. 

Solution  Results.  The  suboptimal  control  program  described  in  the 
preceding  section  was  run  with  the  following  set  of  specified  initial  and 
final  conditions.  The  initial  radius  and  inclination  were: 

r^  =  4263.0  miles 
i^  =  10.0  degrees 

The  final  radius  and  inclination  were: 

r  =  4433.057  miles 
nom 

i _  =  10.746  degrees 

nom 

The  specified  final  radius  and  inclination  are  the  results  of  a  computer 
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program  vhich  calculated  the  Ar  and  the  Ai  of  1000  orbits  with  the 
following  w(r)  control  profile. 

ouo-). 


$•00  tooo 


No.  o  ^ 
OtW\\s 


Fig  8:  Control  Profile  for  Specified  Final  Conditions 

This  control  profile  corresponds  to  the  case  where  the  vehicle  does 

500  radius  change  orbits  and  then  500  inclination  change  orbits.  This 

is  called  the  bang-bang  orbital  transfer  maneuver. 

The  suboptimal  control  program  was  started  at  900  orbits  with  a 

quadratic  control  equation,  './hen  the  program  converged  on  a  set  of 

coefficients  which  matched  the  end  conditions.  r„  equal  to  r  and 

*  x  ^  non 

i„  equal  to  i  ,  and  satisfied  Fq  28,  the  number  of  orbits  was 
i  nom*  ^ 

decreased  by  one  and  the  program  started  again.  The  rationale  behind 
this  strategy  was  that  the  coefficients  of  the  n  orbit  case  would  be 
good  initial  guesses  for  the  (n  -  1)  orbit  case.  It  was  discovered 
that  after  a  few  orbits  had  been  subtracted  from  the  900,  the  program 
refused  to  converge.  At  this  point  another  coefficient  was  added  to  the 
control  equation.  Again  the  program  ran  smoothly  until  a  number  of  orbits 
had  been  subtracted.  Since  the  object  of  this  program  was  to  prove  the 
existence  of  a  solution  to  the  example  problem  stated  in  the  first  section 
of  this  chapter,  the  program  was  terminated  after  the  coefficients  for 
880  orbits  had  converged.  The  converged  coefficients  for  880,  890,  and 
900  orbits  were  chosen  to  show  the  results  of  the  program. 
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The  three  cases,  that  were  chosen,  differed  by  the  number  of  orbits 
and  the  number  of  coefficients  in  the  control  equation.  The  three  cases 
are  summarized  in  Table  1,  The  coefficients  for  each  of  the  three  cases 
can  be  found  in  Table  2.  Figure  9  shows  the  plots  of  v(r)  vs  r  for  the 
above  three  cases.  Table  3  shows  the  savings  in  tine  between  the  three 
cases  chosen  to  represent  the  solution  of  the  total  orbital  transfer 
problem  and  the  nominal  values  for  the  1000  orbits.  These  savings  in 
time  are  equal  to  fuel  s'-ings  since  the  thrust  of  the  vehicle  was 
assumed  to  be  at  full  thrust  throughout  the  entire  problem. 


Case  No.  of  Orbits  No.  of  Coefficients 

1  900  3 

2  890  4 

3  880  7 

Table  1:  Sample  Results 
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Case 

\) 

A1 

A2 

A3 

A4 

A5 

A6 

1 

-.23051 

.24129 

1.41085 

0 

0 

0 

0 

2 

-.05269 

.183179 

1.23758 

-4E"6 

0 

0 

0 

3 

.073167 

.13706 

1.1349 

-.0266 

1 

• 

O 

— w 

O 

_ k 

-.0044 

-.0045 

Table  2:  Control  Equation  Coefficients 


Case 


%  decrease 


1 


10.44 


2  11.53 

3  12.64 


Table  3:  Fuel  Savings 
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V  Conclusions  and  Recommendations 

The  following  conclusions  and  recommendations  are  based  on  the 
results  of  this  study. 

The  equation  derived  for  the  coplanar  circular  radius  change,  Eq  10, 
is  no  more  complex  than  the  equation  given  in  Reference  3  (l69)  for  the 
instantaneous  high  thrust  case.  Equation  10  lends  itself,  very  easily, 
to  computer  use. 

The  equation  derived  for  the  inclination  change  problem,  Eq  39,  is 
also  no  more  complex  than  its  counterpart  in  Reference  3  (169).  This 
equation  also  is  easily  implemented  on  the  computer. 

The  third  conclusion  of  this  study  is  that  significant  savings 
(10%  to  12%)  in  fuel  and  time  can  be  achieved,  over  the  1000  orbit 
bang-bang  maneuver,  if  a  simple  control  profile'  is  followed  in  the  total 
orbital  transfer  problem.  This  savings  in  fuel  can  be  translated  directly 
into  an  increase  in  payload  weight.  Therefore  increases  in  mission 
capabilities  can  be  achieved  without  a  loss  of  vehicle  performance  or  an 
increase  in  mission  time. 

There  are  three  areas  which  can  be  expanded  on  from  this  study.  The 
first  is  to  find  an  analytical  equation  for  w(r)  so  that  the  general  orbit 
transfer  problem  can  be  solved.  The  second  is  the  use  of  a  thrust  angle 
which  allows  the  full  thrust  to  do  both  radius  and  inclination  changes  at 
the  same  time.  This  paper  dealt  only  with  an  angle  that  allowed  either 
radius  or  inclination  changes,  but  never  both.  Probably  the  most  important 
area  in  which  to  investigate  is  the  use  of  complex  orbital  dynamics.  This 
paper  dealt  only  with  two  -  body  dynamics.  Two  -  body  dynamics  is  a  good 
approximation  for  low  earth  orbits,  but  as  the  radius  of  the  orbit  is 
increased,  the  accuracy  of  this  approximation  decreases.  Therefore,  for 
missions  to  geosynchronous  orbit,  complex  orbital  dynamics  must  be  used. 
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